/*
 * Project CHEN Rendering Engine : /optics.hpp
 * Created: 2020-12-16, Last modified:  2021-03-07
 * 
 * There are the most common used optics functions.
 * 
 * Copyright (C) 2021 Hineven (hineven@pku.edu.cn)
 * 
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

#ifndef _H_OPTICS
#define _H_OPTICS
#include "heads.hpp"
#include "geometry.hpp"
namespace chen {

class Fresnel {
public:

	Float etaNegative, etaPositive;
	// 'positive side' is determined by surface normal
	Fresnel(Float erela): etaNegative(erela), etaPositive(1) {}
	Fresnel(Float eneg, Float epos): etaNegative(eneg), etaPositive(epos) {}

	// Return the energy reflection propotion.
	inline Float evaluate(Float cosIn, Float cosTran, bool swap) const {
		Float x = etaPositive, y = etaNegative;
		if(swap) std::swap(x, y);
		Float xa = x*cosIn, xb = x*cosTran;
		Float ya = y*cosIn, yb = y*cosTran;
		Float Rs = (xa-yb)/max(xa+yb, epsilon);
		Float Rp = (ya-xb)/max(xb+ya, epsilon);
		return (Rs*Rs+Rp*Rp)/2;
	}

	inline Vector3 transmitDirection (const Vector3 & w) const {
		Float costheta = w.z;
		Float sintheta = std::sqrt(max((Float)0, (Float)1-costheta*costheta));
		Float singamma;
		if(w.z < 0) singamma = sintheta*etaPositive/etaNegative;
		else singamma = sintheta*etaNegative/etaPositive;
		Float s = pw2(singamma);
		Float cosgamma = std::sqrt(max((Float)1-s, (Float)0));
		Float sqrs = std::sqrt(s);
		Float nm = std::sqrt(max((Float)0, (Float)1-w.z*w.z));
		if(nm != 0) sqrs = sqrs/nm; 
		return Vector3(sqrs*w.x, sqrs*w.y, (w.z < 0) ? -cosgamma : cosgamma);
	}
};

}
#endif
